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Abstract 

We have studied the transition in shape of two dimensional polyampholytes us- 
ing Monte Carlo simulation. We observe that polymers with randomly charged 
monomers get into a globular shape at lower temperatures, provided that their 
total charge is below a critical value Q c . Collapse into globular form happens 
for all forms of force law, but the critical charge depends on the form of the 
force law, inversely dependent on the range of interaction. The value of the 
critical charge is proportional to N, the size of the polymer, in the thermody- 
namic limit. 
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I. Introduction 



The shape of polymers in equilibrium, determines some of its physical, chemical and bio- 
logical properties, thus predicting the shape for a given chemical composition of a polymer is a 
challenging problem of polymer science. Some of the interesting analytical approaches to this 
problem come from the theory of spin glasses 0,0]? the replica trick |^-|5[], the random energy 
model and the mean field approximation p|-[T0[ . 



Due to extremely complicated nature of forces involved, computer simulations have been the 
more successful tool in this study []TT |. Modeling polymers on basis of a self avoiding random 



walk (SAW) using computer simulations can predict some elementary geometrical quantities 
such as the end to end distance (R e ) or the radius of gyration (R g ). For some large polymers 
these two quantities have the same scaling with respect to N, the number of monomers, 

< R 2 e >~< R] >~ N u . (1) 

The universal scaling index u, can be obtained through Monte Carlo simulations, which agrees 
well with Flory's calculation of v for different dimensions, for a neutral SAW JT2|-|il| (table 1). 



Flory's calculations of v are exact expect in d — 3, which never the less is in good agreement 
with simulation (table 2). 

Real polymers however have interactions among their monomers and between the monomers 
and the solvent. Thus they behave differently to the predictions of table 1. However all 
polymers in a good solvent, independently of their monomer interaction, behave like a SAW 
at high temperatures. Thus observed values of v at high temperatures should correspond to 
those of table 1 |T3. In low temperatures the effect of forces become important, thus we expect 



different values of v. 
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In low temperatures the polymer should get into configurations which minimize the potential 
at the same time as satisfying the self avoiding constraint. For short range attractive forces, 
such as the Van der Waals forces, the polymer gets into a d dimensional shape, where d need 
not be an integer. The index v in this case is 1/d. 

This phase transition, which qualitatively happens at the point where the thermal energy 
is comparable with the potential energy, is referred to as the #-point. This is the transition of 
a polymer from a free globular shape to a collapsed globular shape, determined by the minima 
of the potential. In this paper we study the 9 transition for two dimensional polyampholytes, 
polymers with randomly charged monomers. In a random distribution of charge the total charge 
of the polymer need not be zero. If the net total charge becomes proportional to the size of 
polymer, we then expect a behavior similar to polyelectrolytes, polymers with singly charged 
monomers. Higgs and Joany [T9|] have shown that polyampholytes in a good solvent, which 



have total vanishing charge have a smaller radius of gyration than uncharged polymers. They 
have estimated the size of the polymer assuming screening of the long range electrostatic force 
between monomers and using the Debye-Huckel estimation of the free energy. Polyampholytes 
with electrostatic interaction and without the screening effect, were studied by Kantor et al. 
20| , |2l[| , in three spatial dimensions. Using scaling arguments and Monte Carlo simulations, they 



showed that polyampholytes with zero net total charge, collapse into a globular form below a 
critical temperature (T c ). They also showed that a critical net charge Q c exists below which 
the polymer collapses, just like uncharged polyampholytes, where as above Q c the polymer gets 
into a rod shape just like a polyelectrolyte. Of course for T > T c both types, irrespective of 
their total charge, behave like a neutral SAW. 

We have addressed this problem for polyampholytes in two dimensions. This question 
may be relevant when considering polymers in a thin film solvent, adsorbed onto a surface. 
The restriction of movement in two dimensions, do give different results for a SAW in two 



dimensions. But we do not expect qualitatively different behavior with regards to temperature 
or total charge. Indeed the same qualitative behavior is observed, below a critical temperature 
collapse into globular form happens. There also exists a critical total charge Q c , below which, 
the polymer collapses just like a zero net charge polyampholytes, and above which a rod shape 
is taken, similar to polyelectrolytes. 

However the interesting question is the dependence of Q c on the form of the potential. For 
short range potentials Q c is proportional to N, the size of the polymer. We performed our 
simulations for a few different forms of the potential, and find agreement with an approximate 
expression for Q c 

Ql ~ {aN 2 + bN)q 2 . (2) 

Where q is the charge of the monomer, iV is the size of the polymer, a and b are potential 
dependent constants. We give an approximate derivation of eq. (2), which is not dependent on 
the dimension of space, so it should also hold in higher spatial dimensions. 

II. Simulations 



To simulate the motion of a polymer we have adapted Rose's model ||22|| . In this model 
a polymer, made up of N monomers is replaced with a SAW with N — 1 steps on a square 
lattice. Every vertex of the lattice through which the walk has passed, represents a monomer. 
Since two monomers can not occupy the same location in space, a SAW is a representation of a 
configuration of the polymer. To allow for effects of interaction, each site of the lattice on the 
random walk is allowed to take a random charge, which we restrict to be ±q for simplicity. 

To calculate quantities such as the radius of gyration one needs to perform averages over 



ensembles of configurations. We used the Metropolis JTTJ] method of importance sampling to 
generate these ensembles. Let us assume C is the set of all possible configurations. We then 
construct a set of maps f a which transform the members of C into each other: 



f a : C - C 



fa [ s i 



Si, Sj G C. 



(3) 



At any given moment, a weight factor determines the probability of the action of f a , using the 
interaction hamiltonian. Thus starting with a given configuration, the action of f a takes us 
through the phase space C with the passage of time, generating the more probable configura- 
tions first. An interesting question is what is the most efficient set of f a , such that any two 
members of C can be transformed to each other with a finite number of actions of f a . A very 
large set of f a is not desirable since the computer has to run through the set at each time step, 
which slows down the algorithm considerably [22]. On the other hand a very small set of f a 
leaves some configurations unreachable. Consider the following set of moves: 

i - Standard move 



ii - End moves 



(a) 



(b) 
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Hi - Crank shaft move 



• • 



• • • • 



• • • • 



• • 



iv - End reflection 



• • 



• • • • 



• • • • 



• • • 



Without the moves (ii.b) and (iv), the following configuration: 

• • • 



can not be reached. However, the other moves may be sufficient for simulation of SAW, at high 
temperatures because the compact configurations such as the above have very small entropy. 
But for our problem, where we seek compact configurations, at low temperature we need the 
full set of above moves. Apparently this set of moves seems to be sufficient for our problem, 
however we don't have a proof of this. 

For a given polymer of size N, the algorithm running time grows exponentially with N. For 
each configuration iV attempts at changing the configuration is made. The number of possible 
configurations grows like N 2 in two dimensions. Note that the algorithm has to check the illegal 
moves as well. Thus the time to reach equilibrium grows like N 3 . In practice we observed that 
approximately 100A r3 was necessary to reach equilibrium. On the other hand, calculating the 
change in energy, going from one configuration to another, requires N calculations due to the 
long range nature of the force law. Thus run time grows like iV 4 , hence very long polymers are 
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not easy to simulate. We were thus restricted to a size N < 32, with once exception where we 
succeeded in getting result for iV = 48, too. 

Attributing a fixed charge ±g to monomers at random, leaves us with a polymer with a total 
random charge, between +Nq , and —Nq . We then choose a monomer at random and change 
its charge until the desired total charge is achieved. This gives the initial configuration to start 
with. We then allow this configuration reach equilibrium for a given temperature following the 
Metropolis method [jn|| , applying the moves already discussed. After reaching equilibrium the 
radius of gyration is calculated, and this is then repeated 30 times at large time intervals. The 
whole procedure is then repeated with different charge configurations. Averaging over charge 
configuration a final value for < R g > is obtained. Repeating the above for different value of 
N gives a scaling relationship between < R g > and N, giving a value for v. 

III. Results 

We simulated polyampholytes with size N = {32, 16, 8, 4} in two spatial dimensions, using 
different force laws acting between monomers. We assumed an electrostatic interaction of form 
log(r), -3j 3 and a short range interaction U(r) = for r > I, and constant otherwise. To 
begin with we simulated polymers with total charge zero, and observed that there exist a 
phase transition as in three dimensions, where the polymer gets into a globular shape at low 
temperatures. This phase transition is observed for all types of potentials used. 

In figs, (l.a) and (l.b) the radius of gyration is plotted vs. [3 = (™), for the two potentials 
log(r) and -. Where we observe a reduction in R g as (3 increases. This reduction is better seen 
for large polymers. We can observe the changes in v as a function of f3 for these two potentials in 
tables (3. a) and (3.b) and figs. (2. a) and (2.b). We observe that this two types of polymers have 
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similar behavior and the value of u for high temperatures is consistent with the two dimensional 
neutral SAW. The interesting point is that lower temperatures, these polymers cannot collapse 
completely, thus v > 0.5, which is the index for a two dimensional shape. This may have been 
caused by considering an incomplete set of moves, which can not access a collapsed shape. 
However increasing the number of moves and repeating the simulations and consistency of our 
results suggest otherwise. In two dimensions holes can form which ore surrounded by like 
charges, thus, movement to a more compact configuration becomes energetically inaccessible. 



Probability of having a length I of like charged monomers is proportional to a power of / [p3 
thus the distribution and size of these holes obeys a scaling law. Therefore the existence of such 
holes affects the scaling dimension of globular collapse in two dimensions and forces v above 
0.5. With purely attractive forces these holes would not exist and this obstruction is removed. 

Repeating simulations for polymers with total nonzero charge, showed that the phase tran- 
sition observed by Kantor and Kardar [BT| in three dimensions, exists in two dimensions as 
well. For a total charge greater than Q c , rod shape is observed where as for Q < Q c , globular 
form, just like a zero net charge is observed. The interesting point is the dependence of Q c on 
the potential and the size of the polymer. 

In figs. (3a-3c) we observe the results for < R g >, for a polymer of length N = 32, using 
different potentials, and in figs. (4a-4c) results for N = 48, 32, 16, using the potential of the eq. 
(12). Although a critical Q c dose exist for each potential, but the value of Q c depends inversely 
on the range of the potential. 

IV. The Critical Charge 



The charge on monomers is a random variable of the problem but without dynamics, thus 
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it should be treated as a quenched variable when analyzing this systems. This is the root of 
similarity between this problem and spin glass systems. The appropriate quenched average is 
to be effected as the free energy 



F = log Z (4) 

where bar means average over quenched variable, in our case the different charge configuration. 
The calculation of quantities such as eq. (4) is riddled with difficulties, one way is the replica 
method 0fj], but there is a lot of controversy associated with it. We were unable to perform 
a direct calculation of eq. (4) for our system. Let us instead often a compromise quantity the 
weighted average potential energy: 

« U »= J J D r P(r)D q P(q)Y^U(r t ,r 3 )e~P u M. (5) 

hi 

In this expression D r P(r) integrates over all possible self avoiding configurations including 
the appropriate distribution in space, D q P(q) like wise averages over all charge configurations, 
the factor e~@ u allows the more probable energy configurations to happen more frequently in 
the averaging. The expression of eq. (5) lacks a normalizing factor which can be added later. 
The sum over i and j, sums over all monomer pairs in interaction. Taking out the weight 
factor e~^ u we end up with the expression used by Kantor and Kardar [21|. In this case the 
calculation of eq. (5) leads to: 

« U »~ (Q 2 - Nq 2 Q ) (6) 

which indicates that the average potential changes sign at the point of total charge Q equaling 
q y/~N. Thus a phase transition is expected at critical total charge of q y/N. This was observed 
by Kantor and Kardar but is inconsistent with our observations. 
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Not including the Boltzman factor implies that two spatially similar configuration, with 
different charge configurations are equally probable. Where as clearly when opposite charges 
are near each other, a greater probability is achieved due to lower energy. Let us now keep the 
Boltzman factor in eq. (5). We can write eq. (5) as : 

« U »= n + J D r P(r)U(r)e- mr) - n_ f D r P(r)u(r)e- pui - r) (7) 

where n + and n_ are the number of pairs of monomers with like and unlike charges respectively: 

N(N - 2)ql + Q 2 



n+ = 



4 

N 2 g 2 - Q 2 



n- = ^ (8) 



which results in: 



, 2 2N + (X-1)N 2 2 . 



« U »~ (Q 2 A V +1 ' q 2 ) (9) 



where 



JD r P(r)e-^r { } 

The parameter A, can only be determined if the shape of the potential is given, but A is positive 
and larger than 1. We now observe that the average potential energy changes sign at the critical 
charge: 

Ql = [( XTT )iv2 + xTi N]q "- 

Thus the critical charge is proportional to N, rather than y/N. This is in better agreement 
with our results. We can calculate A for a simple potential: 

V r<l 

U(r) = { ' (12) 

r> I 
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We have A = e 2l3V . Comparison with fig. (4.b) indicates 16 < Q c < 20 for N = 32, which is 
consistent with eq. (11). This result indicates that in the thermodynamic limit where N — > oo, 
critical charge scales with N. If the factor A calculates to be almost near 1, then we get a 
different scaling behavior where Q c scales like N% . This is of course true for high temperatures, 
i.e. A = 1, but in this region, interactions don't play a role and all polymers behave like neutral 
self avoiding random walks. 

V. Discussion 



Simulating polyampholytes in two dimensions we observed the qualitative feature that a 
globular form is taken at low temperatures if the total charge of the polymer is less than a 
critical value Q c . The scaling index u, of this globular form is above 1/2, this is conjectured 
to be due to formation of areas surrounded by parts of the polymer made up of monomers of 
like charge. Where such an area is formed other parts of polymer are forbidden to enter it, 
due to the self avoiding nature of this system. Thus the globular form is more swollen than 
a two dimensional shape, hence v > 1/2. We were unable to predict the deviation from two 
dimensions by an analytic calculation, although this seems possible. Clearly a potential which 
is only attractive such as the Van der Waals potential would not give rise to this effect. 

The other interesting observation concerns the scaling of Q c with size. Our analysis indi- 
cates that Q c scales with N. This becomes more stark when the forces are short range. It 
would also be experimentally observable if the effects of screening are made more intense, for 



example by changing the solvent p^ , pq| . It appears that this effect should exist independent of 
dimension, although our simulations were restricted to two dimensions. For a long range force 
the dependence of Q c on N may become smaller, but should be important for large enough N. 
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It must be said that our calculations were not based on a true quenched average. Since per- 
forming a true quenched average seems impossible, therefore some more convincing calculation 
is necessary before the question may be settled. 

Acknowledgments: We are indebted to Mehran Kardar for valuable discussions and sug- 
gesting this problem to us. 
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Tables 



Monte Carlo Simulations for d = 3 
0.59992 ± 0.002 g| 
0.5745 db 0.008 ± 0.0056 @ 
0.588 ± 0.001 



table 1. Flory's calculation of v for different dimensions. All are exact except d = 3. 
table 2. Monte Carlo simulation estimates of v for d = 3. 
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0.02 
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0.009 




0.02 
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0.008 


0.08 


0.753 


0.009 




0.08 


0.761 


0.009 


0.32 


0.726 


0.010 




0.32 


0.744 


0.010 


1.28 


0.644 


0.013 




1.28 


0.672 


0.017 


5.12 


0.614 


0.016 




5.12 


0.590 


0.032 


20.48 


0.60 


0.012 




20.48 


0.61 


0.036 



table (3.a) table (3.b). 



d 


Flory's v 


1 


1 


2 


3/4 


3 


3/5 


> 4 


1/2 



Table 3. The dependence of v on (3 for polymers with zero net charge for the potentials log(r) 
.a), and i (3.b). 
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Figure Captions 



Figure 1. 

The radius of gyration of the polyampholyes with zero net charge and 
different numbers of monomers vs. (3 = l/T for the potentials log(r) 
(l.a), and ± (Lb). 

Figure 2. 

The dependence of v on (3 for polymers with zero net charge for the 
potentials log(r) (3. a), and £ (3.b). 

Figure 3. 

The radius of gyration of the polyampholyes with net charge Q and 
N = 32 vs. (3 for the potentials log(r) (3. a), £ (3.b), and ^ (3.c). 

Figure 4. 

The radius of gyration of the polyampholyes with net charge Q for the 
potential of eq. 12 vs. (3 for TV = 16 (4.a), N = 32 (4.b), and TV = 48 
(4.c). 
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